use vars qw($USAGE);

use strict;
use Getopt::Long;
use Bio::SeqIO;

$USAGE = "longorf [--help] [--notstrict] [--verbose] [--graph] [--width printwidth] [--format seqformat] --input seqfile\n";

my ($sequencefile,$sequenceformat,$notstrict,$graph,$verb,$printwidth,$help) =
   (undef, 'fasta', undef, undef,undef,50, undef);

&GetOptions('input|i=s' => \$sequencefile,
            'format|f=s' => \$sequenceformat,
            'notstrict|n' => \$notstrict,
            'width|w=s' => \$printwidth,
            'graph|g' => \$graph,
            'verbose|v' => \$verb,
            'help|h' => \$help,
            );

if ($help) {
   exec('perldoc', $0);
   die;
}

if (!defined $sequencefile) {
    die($USAGE . "\nPlease specify an input filename.\n");
}

sub longestORF {
   my $best=0;
   my ($bests,$beste,$beststrand)=(-1,-1,0);
   my $bestorf="";

   my $relaxed=$_[1];
   my $dna=Bio::Seq->new(-seq => $_[0]);
   #my %strand=('+'=>$dna->seq,
   #            '-'=>$dna->revcom->seq);
   my %strand=('+'=>$dna->seq);

   foreach my $direction (keys %strand) {
      my @starts=();
      my @ends=();
      if ($relaxed) {
         for (my $frame=0;$frame<3;$frame++) {
            unless ($strand{$direction}=~m/^.{$frame}(taa|tga|tag)/i) {
               push @starts,$frame+1;
            }
         }
      }
      while ($strand{$direction}=~m/(atg)/gi) {
         push @starts,pos($strand{$direction})-2;
      }

      while ($strand{$direction}=~m/(taa|tga|tag)/gi) {
         push @ends,pos($strand{$direction})-2;
      }
      push @ends,($dna->length-2,$dna->length-1,$dna->length);

      for my $s (@starts) {
         for my $e (@ends) {
            if ($e%3==$s%3 and $e>$s) {
               if ($e-$s>$best) {
                  $best=$e-$s;
                  ($bests,$beste,$beststrand)=($s,$e,$direction);
                  $bestorf=Bio::Seq->new(-seq=>$strand{$direction})->subseq($s,$e);
               }
               last
            } else {
               next
            }
         }
      }
   }
   return ($best,$bests,$beste,$beststrand,$bestorf);
}

my $seqio = new Bio::SeqIO('-format' => $sequenceformat,
                           '-file' => $sequencefile );

my ($length,$start,$end,$direction,$sequence);
my $count=0;
my @lengths;
my $totallength=0;

while (my $dna = $seqio->next_seq) {
   $count++;
   ($length,$start,$end,$direction,$sequence)=longestORF($dna->seq,$notstrict);
   if ($verb) {
      print $dna->display_id," ",$dna->desc,": ";
      print "$length, $start, $end ($direction)\n$sequence\n\n",Bio::Seq->new(-seq=>$sequence)->translate->seq,"\n\n--\n\n";
   }
   $totallength+=$length;
   $lengths[$length/3]++;

   # print the longest orf length for each sequences
   print $dna->display_id, "\t", $length, "\n";
}